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FOREWORD 


The work described in this report was performed for the Auxiliary 
Propulsion Branch of the Propulsion Division, Structures and Propulsion 
Laboratory, NASA George Q. Marshall Space Flight Center under Contract 
NAS8-33974t The NASA technical monitors were Lee W. Jones and T.D* McCay. 

This report represents Part II of the Final Report. Part I deals with 
the evaluation of alternate laser devices as well as alternate propellants 
and energy coupling schemes. 


ABSTRACT 


The present status of theoretical models for the laser heated thruster 
is reviewed. It is concluded that existing models neither agree with each 
other nor with the limited experimental data available. The requirements 
for an improved laser heated thruster theoretical model are discussed. The 
application of a time-dependent finite-difference Navier-Stokes equation 
solution to the laser heated thruster problem is described, along with a 
simple closet^, form solution which was developed in order to gain insight 
into the difficulties encountered in the pursuit of the numerical solution. 
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1. LASER HEATED PROPULSION MECHANISM 


The concept of laser heated propulsion is based on using a laser beam 
to remotely heat a working medium in a rocket thruster. Such a scheme 
appears to be attractive because it might be possible to generate propellant 
gas temperatures which are much higher than those possible in chemical pro- 
pulsion devices, in particular, in conjunction with propellants of low 
molecular weight, a very high thrust chamber temperature will result in a 
specific impulse much higher than that achievable in combustion driven 
rocket motors. Even through the concept is simple, the attendant high 
temperature presents difficult engineering problems. 

The use of continuous wave laiiers for laser heated propulsion can be 
divided into a very high specific impulse, high temperature regime and a 
high specific impulse, medium temperature regime. While the first regime 
would be realized by utilizing hydrogen as the propellant and a laser 
supported absorption wave generated and sustained by inverse bremsstrahlung 
absorption, in the medivtm temperature regime laser energy would be coupled 
into a suitable propellant by molecular or aerosol absorption. Details of 
the'^e and other schemes are discussed in Part X of this report. 

The work described here is mainly coricerned with an assessment of laser 
propulsion fluid mechanics modeling and simulation requirements for the very 
high specific impulse, high temperature regime. Here the essential feature 
is that the propellant gas, which is normally a transparent, non-conducting 
medium, becomes ionized by the laser beam and then strongly absorbs energy 
from the beam (Ref. 1). 

If a few free electrons are present in the focal volume, they can gain 
sufficient ene)rgy from the beam to produce further ionization via collisions 


with neutral gas particles. Breakdown will occur if the electron density 
reaches a critical value despite losses due to diffusion, attachment, 
recombination and other mechanisms. 

Once a plasma has been formed, it emits radiation whose wavelength dis- 
tribution depends on the temperature and the density. Host of the energy 
emitted is in the far ultraviolet and is believed to be more readily 
absorbed in the surrounding cold gas rather than in the plasma itself. 

Thus, a layer of gas outside the plusma, although transparent to laser 
radiation, is heated by the plasma radiation, and upon reaching a suffi- 
ciently high temperature, will be ionized to such an extent that it will 
also become absorbing for the laser beati>. This layer will then be further 
heated rapidly until its temperature becomes so high that it in turn will 
become transparent to ultraviolet plasma radiation. By this time a new 
layer of plasma nearer the laser will have become absorbing, so the boundary 
of the plasma, called an absorption wave, will move toward the laser. If 
this occurs in a flow of propellant gas, the flow velocity of which is equal 
^ III magnitude but opposite in direction to the motion of the absorption wave, 
we have a stationary propellant plasma' which can function as the heating 
mechanism for a propulsion device. 

Various laser heated propulsion devices have been conceived and 
proposed. Since it would be difficult to discuss the modeling and simula- 
tion requirements in complete generality, one particular configuration for a 
laser heated thruster is shown in Fig. 1 (schematically). This is the 
arrangement used in the subsequent discussion of modeling and simulation 
requirements. In particular, we are interested in describing or modeling 
the interaction of the initially cold propellant gas (H2) as it flows 
through and around the laser supported plasma region, then mixes and exits 
through the subsonic-supersonic nozzle, producing thrust in the process. 
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^ 1 - Laser Heated Propulsion Device (Schematic) 





Fur tha actual. Uesigu of 8uch a device it is most important to predict 
the propellant velocity required to achieve a stable plasma as well as the 
cooling requirements for the plasma and mixing chamber walls. Hydrogen is 
of particular interest as a propellant because it promises the highest 
specific impulse. 
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2. REVIEW OF EXISTING MODELS 


2.1 ABSORPTION WAVE MODEL | 

Although much has been achieved in modeling the physics o£ the laser | 

heated thruster, much must be done toward obtaining agreement with experi- | 

ment if these mathematical models are to serve as useful tools in the design I 

of laser heated propulsion devices. Most of the existing models deal with | 

absorption waves in air and have numerous simplifying assumptions in common. | 

Very little work so far has been done with hydrogen as the propellant. j| 

Experimental data for absorption waves in hydrogen are also very sparse. | 

fi 

The first model to appear in the literature was that of Raizer (Ref. f 

2). His model is based on the assumptions of one-dimensional flow at 
constant pressure with heat conduction concidered to be the principal heat 
transfer mechanism. The plasma is considered to be optically thin to its ; 

self-radiance, and in thermal equilibrium. This allows the electron density 
to be determined from the Saha equation. Also, the absorption coefficient 
for the absorption of laser radiation via inverse bremsstrahlung as well as 
an approximate function for the overall radiation losses can be specified as 
functions of temperature and pressure. 

Finally, in order to arrive at an analytical solution, the absorption : 

coefficient as well as the ratio of specific heat to thermal conduction are | 

assumed constant. The neglection of radiation losses limits the applic- | 

ability of the solution to small plasma volumes (diameters less than a few 
millimeters at p » 1 atm). 

With the goal to improve upon the analysis by Raizer, Jackson and 
Nielsen (Ref. 3) presented a model which includes the mechanism of radiative 
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iransjEer into the orobient gas (air)* Apart from this, other assumptions 
made are basically the same as those invoked by Raiser* Jackson and Nielsen 
derive an analytic solution, the evaluation of which is performed by iterat" 
ing on the value of the radiative, loss function. The authors point out 
several shortcomings of their analytical solution: (1) interferograms taken 

of absorption waves show clear evidence of radial flow, rendering the prob- 
lem two-dimensional, and (2) their calculated propagation velocity of the 
absorption wave as a function of incident laser intensity does not agree 
with experiment. 

Shortly following the analysis of Jackson and Nielsen, Batteh and 
Keefer (Ref. 4) presented a two-dimensional generalization of Raizer's 
analysis for the subsonic propagation of laser sparks in air. To be pre- 
cise, Batteh and Keefer combined a one-dimensional flow velocity description 
with a two-dimensional temperature field description, retaining Raizer's 
assumption of ordinary thermal conduction as the principal propagation 
mechanism, with plasma radiation serving only as an energy loss mechanism. 
The importance of this two-dimensional temperature field model is that it 
provides insight into the effect of boundaries on spark propagation in a 
channel. Other assumptions such as constant pressure, constant ratio of 
specific heat to thermal conductivity are also retained. The absorption 
coefficient was assumed to be piecewise constant, zero in front of the 
absorption wave, and non-zero within the plasma region. An analytic solu- 
tion in terms of Bessel functions was obtained to evaluate the relationship 
between laser intensity and spark propagation velocity. 

Batteh and Keefer found that their calculated propagation speeds are 
roughly an order of magnitude lower than those obtained experimentally, 
while the maximum temperature in the wave as calculated matched that of the 
experiment. The authors suggest that this is due to considering thermal 
conduction only as the primary propagation mechanism. They essentially 
argue that plasma radiation, causing ionization and additional laser energy 
absorption via inverse bremsstrahlung in front of the absorption wave must 
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be included in the model to achieve realietic prediction! of the propagation 
velocity ai a function of laaer intenaity* 

The next model, that of Kemp and goot (Kef. 5), in fact doea include 
the effect of forward plaama radiation in a one>^ime(iaional numerical 
model. Thia model aaaumea one-dimeriaional flow at conatant preaaure in a 
constant area channel. Since the equations are solved by numerical integre- 
tion, it is possible to treat the thermodynamic properties as well as the 
various transport phenomena (absorption, thermal conduction and radiation) 
in a much more detailed manner than in the previously discussed analytical 
solutions. The solution of Kemp and Root's formulation requires an itera<~ 
tion on the mass flux for specified laser intensity. Only one particular 
mass flux value will yield a ptysically reasonable temperature distribution 
in the absorption wave (saddle-point singularity). Calculations are pre- 
sented for pressures of 1, 3, 10 and 30 atm, and for laser powers of 10 kW 
and 3 MW. The authors observe a very steep rise in temperature at the 
leading edge of the absorption wave and caution that for such conditions the 
"radiation conduction approximation," used to consider the heating effect of 
plasma radiation on the cold gas in front of the wave, is not really valid. 
However, since energy is properly conserved in the model, they suggest that 
Che temperature profile in the wave proper should be qualitatively correct. 
No experimental data for absorption waves in hydrogen were available at the 
time Co confirm the theoretical predictions. Also, the authors realize that 
the one-dimensional nature of the model presents serious limitations on the 
validity of the results. It is pointed out Chat thermal conduction in the 
radial direction becomes important for small plasma regions at low pressure. 
A very obvious two-dimensional phenomenon is Che focused laser beam. In 
this case, the intensity will not only change due to absorption but also due 
to varying cross-sectional area. A convergent beam behaves as a higher 
intensity beam - relative to a collimated beam of the same initial intensity 
- and requires a higher propellant flow velocity to keep the "plasmatron" 
sCaCionarr. Finally, there is Che transverse velocity component. It is 
suggested that its effects should be less in a confined flow - such as in a 
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thrust^er - than in thn unconfined caae. Since neither theoretical nor ex- 
perimental reaulta are availablei it ia difficult to aaaeas the importance 
of two-dimenaional flow for abaorption wavea in a channel. 

Xn concluding thia brief diacuaaion, it appeara that the modeling of 
abaorption wavea ia atill in a atage of infancy. Xn Fig* 2 (adopted from 
Ref. 3) we have aummarlzed reaulta for abaorption wavea in air from the work 
reviewed here. Some experimental reaulta are alao ahown. No attempt haa 
been made to aystematically correlate theae reaulta. It ia evident, how- 
ever, that theoretical predictiona are in conaiderable diaagreement with the 
experimental data (Refa. 6^7) and with each other, even for aimilar aaaump- 
tiona in the theoretical modela. 

2,2 THRUSTER AND NOZZLE FLOW MODELS 

Xn order to produce thruat, the heated gaa exiting from the abaorption 
wave inuat be accelerated through a auperaonic nozzle which baaically 
converts thermal energy into kinetic energy. Just as for the abaorption 
wave, the modeling of this flow can be approached at different levels of 
sophistication. 

The first of such models to be discussed here ia that of Kemp and Root 
(Ref. 5). It is a one-dimensional model, obtained by explicitly including 
the area term and retaining the momentum equation in the set of one- 
dimensional equations used for the one-dimensional absorption wave model by 
the same authors. The results of the absorption wave calculations serve as 
initial values for the thruster problem, i.e., the solutions are joined at 
an appropriate location downstream of the wave. The model described in Ref. 
5 represents an •'inverse method*" In Ref, 5, rather than specifying the 
thruster and nozzle contour, the velocity is specified as a function of flow 
distance. While this or variations of this procedure are common in order to 
avoid the sonic singularity at the throat, it means that many iterations may 
have to be performed until a reasonable or specified nozzle contour is 
obtained. 


8 



Ref. 4 (Batteh and Keefer) 
/ (p “ 1 a^m) 




While the one-dimensional thruster flow model has the advantage o£ 
simplicityi it cannot properly account for lateral heat condution and radia- 
tion to the walls of the thrust chamber. As pointed out in Ref, 5, a pre- 
cise evaluation of the transverse heat loss due to thermal conduction and 
radiation requires a two-dimensional model. Such a model was presented by 
the same author (et al.) in Ref. 8. This model is based on the boundary 
layer equations and thus includes radial gradients. Added terms for laser 
energy absorption, equilibrium chemistry and plasma radiation losses 
complete the equations. Although this model is an improvement over the 
previously mentioned one-dimensional model , it does not contain any terms 
for forward energy transport by radiation or conduction, and therefore 
cannot be applied to model the absorption wave. A method to join the one- 
dimensional absorption wave calculations to the two-dimensional thruster 
calculations is described in Ref. 9. According to the authors, the two- 
dimensional thruster code needs a number of refinements to make it generally 
useful. Required are a generalization of the radiation model, more exper- 
ience or a rationale to choose the axial pressure gradient to obtain a 
desired flow channel contour, and a variable coordinate grid structure to 
concentrate grid points in areas of steep gradients. 

2.3 THERMODYNAMIC AND TRANSPORT PROPERTIES 

As compared with other fluid flow problems, the description and 
analysis of the flow in a laser-heated thruster is severely complicated by 
the large temperature ranges over which thermodynamic and transport phe- 
nomena have to be considered. As a consequence, the early investigations, 
in particular those that produced analytical solutions, had to employ many 
simplifying assumptions (Refs. 2, 3 and 4)> While the qualitative nature of 
the results from these analyses is useful in furthering an understanding of 
the phenomena involved, the pursuit of quantitative results probably demands 
numerical solutions and much more detailed thermodynamic and transport prop- 
erties, such as used by Kemp et al. (Refs. 5, 8 and 9). In fact, a substan- 
tial portion of the work by Kemp et al. deals with the development, evalu- 
ation and assembly of supporting data such as laser absorption coefficients, 


10 


as v;ell as thermodynamic and transport properties (conduction and radia- 
tion) . It is anticipated that any future modeling effort will rely heavily 
on that work. 

Vflille all the existing models so far have assumed thermodynamic equi- 
librium conditions, it might be necessary to consider nonequilibrium condi- 
tions in the future, as pointed out by Batteh and Keefer (Ref. 4) in their 
discussion of the discrepancy between observed and predicted absorption wave 
propagation velocities. A calculation of radiative properties of nonequi'i- 
librium hydrogen plasma has been presented by Park (Rsf. 10) in terms of a 
computer code which calculates emission and absorption coefficients from 
given electron temperature, electron density, neutral particle density and 
intensifies of incident radiation. In calculating radiative transport the 
program shows that there is a large difference between calculated intens- 
ities of radiation emitted by a bulk of equilibrium and nonequilibrium 
hydrogen plasma. 
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3. THRUSTER MODEL REQUIREMENTS 




Based on the preceding discussion of existing models it is now possible 
to establish the basic features which an improved laser-heated thruster 
model should incorporate • These features are: 

I, Axisymmetric , two-dimensional flow so that effects of radial 
velocity components can be assessed 

2t Axial as well as radial heat transfer, both by thermal 
conduction and radiation, and 

3. Converging (or focused) laser beam configuration to 
determine the required mass flow more accurately. 

Requirement (2) has implications which need further elaborations. It will 
be shown that this leads to additional requirements. 

As previously discussed, the propagation of plasma radiation (VUV) 
upstream to the front of the absorption wave will likely lead to a non- 
equilibrium ionization situation. The reason is that electrons have the 
special property that their particle mass is much less that of any of the 
other constituents. Because of the inefficiency of energy exchange in 
elastic ''.ollisions between particles of disparate mass, the electron temper- 
ature may differ appreciably from the temperature of the heavy particle 
species (Ref. 11). Additional equations must therefore be included in the 
model to account for electron mass and energy conservation. 

A second, and probably more serious implication of requirement (2) is 
the fact that the inclusion of upstream influences renders the mathematical 
problem elliptic. This means that a solution cannot be attempted using a 
single pass forward marching scheme. Instead, the problem has to be solved 
considering the entire domain in a time-dependent fashion. 
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The stated requirements render the potential model more difficult to 
evaluate than previous models, particularly so if it is visualized that the 
basic set of equations governing the fluid mechanics (conservation of mass 
and energy plus the Navier-Stokes equations for conservation of momentum) 
will have to be solved numerically. 
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A. MODEL FORMULATION) ANALYSIS AND RESULTS 
A.l DETAILED MODEL FORMULATION 

A major goal of this effort was to investigate the effect of two- 
dimensional fluid mechanics on the conditions in a laser heated thruster by 
utilizing the available Lockheed GIM code (Ref. 12) which numerically solves 
the full Navier-Stokes equations (in conjunction with mass and energy con- 
servation equations). Using this code we can therefore satisfy requirement 
1, and, as far as thermal conduction is concerned, also requirement 2, as 
stated in the previous section. 

The basic equations for axisymmetric flow as used in the GIM code are: 


at ax r ar ^ 


+ G » 0 


( 1 ) 


where U, E, F and G each represent four components for the equations of con- 
versation of mass, momentum in x- and r-directions and energy as follows: 
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The above equations represent four equations for the four unknowns - p, u, v 
and T - if we express the pressure through the equation of state as a func- 
tion of p and T, l>e., 


P « PRT (12) 

which is valid for a perfect gas. The terra Q in the energy equation (see 
Eq.(5)) denotes the net rate of energy addition to the flow which is repre- 
sented by the difference between laser energy absorbed and plasraa radiation 
emitted. Since a primary goal was to compare our results with those 
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obtained from the Batteh and Keefer analyaia (Ref. 4)> we u9e their 
formulation to express Q as 

Q - (Xl^ m k (T « T.) (13) 

where represents the incident laser energy flux of a collimated beam 
and <X is the absorption coefficient. The second term represents a rather 
crude approximation for the energy loss due to plasma radiation, with the 
form of this term chosen so as to make the Batteh-Keefer anaytical solution 
possible (op. cit.). 

At the present time, the GIM code is formulated in terms of constant 
thermodynamic and transport properties. The rather wide range of tempera- 
tures of interest in the laser propulsion problem therefore requires us to 
select suitably averaged quantities for the specific heat, the thermal 
conductivity and the viscosity. In particular, since we want to compare 
results with those of Batteh and Keefer, we want to use the same average 
value for the ratio of specific heat to thermal conductivity that was used 
in their calculations. For consistency, the same averaging procedure is 
applied to the gas constant. 

Figures 3 through 5 show the thermal conductivity (Ref. 5), the 
specific heat at constant pressure (Ref. 13) and the molecular weight (Ref. 
13), respectively, for equilibrium hydrogen at p > I atm as a function of 
temperature. Integrated averages for c^, k and the gas constant R are 
also shown. From these, we can obtain, for use in the GIM code calcula- 
tions, average values for the ratio of specific heats and the viscosity by 
applying standard relationships such as 
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Fig. 5 - Molecular Weight and Gas Constant of Equilibrium Hydrogen 
at p “ 1 atm (Ref, 13) 
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with the letter taken from Svehla (Ref, 14). Finally, the conatant m in the 
radiation lots function has been evaluated by Keefer et al., to be 

m " 14*3 cm"^ (at p ■■ i atm). 

Before we attempt a numerical solution of the full set of equations 
given in this section, it is useful to investigate possible simplifications. 
As it turns out, for the contemplated conditions such simplifications are 
possible. With additional assumptions, this leads to an explicit analytical 
solutior^ which can easily be evaluated. This will be discussed in the 
following section. 

4.2 A SIMPLE ANALYTICAL SOLUTION 


While it is our primary goal to obtain a numerical solution of the full 
set of equations presented in the previous section, obtaining such a solu- 
tion will be facilitated by some advance knowledge of what that solution 
might look like. With this in mind, let us non-dimensionalize the equations 
by introducing suitable reference quantities for all relevant variables. 
Defining dimensionless variables as 


P' = P/Pq 

u', v' = 1,V v/u^ 

. X' = /i, X/Mq 

= X, y/i, t' =: t 

T' = T/T , = g/c T 

c' = c /c 

V v' VO 

k' = k/k^ 


(16) 


and introducing these into the equations, the equations (with primes dropped 
for convenience) take on the following form: 
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Pue to the non-dimens ionallzing process all flow variable groups are 
now of order unity, And we must investigate the order of magnitude of the 
dimensionless parameters Re, Pr and E. While Pr is generally of order unity 
for a‘i , vases, it can be shown that, for the conditions in the laser heated 
thruster, the Reynolds number is also of order unity, while 2 is several 
orders of magnitude smaller. We can conclude that all terms are of equal 
importance in the continuity and momentum equations, but in the energy 
equation, however, the pressure terms, all carrying the E as a multiplier, 
will be orders of magnitude smaller than the convection and conduction terms. 


In pursuit of a simple analytic solution we iwv examine steady state, 
inviscid one-dimensional flow at constant pressure. With these assumptions, 
the governing equations reduce to 
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Recognising khet + p/p ) repreeente the total enthalpy which » for out 
case, conaiatf of mostly thermal energy | we ignore the kinetic energy and 
write 


^ + £ a H<wh » cp T + constant 


( 23 ) 


Knowing that the GIM code presently requires c^ and k to be constant p the 
energy equation reduces to 


^ ^ Q 
dx^ " ■~TT dx “ - Tc 


(2/f) 


where, from Eq* (13) 


Q« al 0““^ - mk (T-T.) 
o 1. 


(25) 


[ 


Of coursep this ordinary differential equation for T(x) is almost the same 
equation used by Raiser (op, cit.) ei ept that we further simplified it by 
introducing an integrated average k and by nf/glecting radial conduction. As 
a matter of fact, we can easily include raduvl conduction in the same manner 
as Raiser did by modifying the constant m in the plasma radiation loss term. 


We can rewrite the differential equation as 


- a 


dT 

dx 


- mT = - 


P «'Otx ^ rr 

T ® -“^l 


and arbitrarily specify that the laser power density P - obeys 

P s 0 for X < 0 
P = P^ for X > 0 


(26) 
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such that X ■■ 0 designates the location of plasma ignition (onset of ioniza** 
tion and laser energy absorption). Solutions are then obtained separately 
for X < 0 and X > 0, and the solutions and T 2 joined by requiring 
that both T and its derivative match at x ■■ 0. Thus, the boundary condi- 
tions take the following form: 


X s 0: 


■^1 = ^2 = To 

dT, dT 

I St 

ax dx 
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- 00 
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+ 00 


1 ' 
dT2 
dx 
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s 0 


Standard procedures for nonhoinogeneous linear differential equations then 
yield for x < 0: 


and for x > 0 : 


where 
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Imposition of the boundary conditions on l 2 (x) generates an additional 
relationship for the laser power density as a function of a, b and of» i.e.^ 


- 2k [m.a(gf a)] b 
b - (a f 2 a) 


(T^-T.) 


(32) 


The above equations show that the solution has one free parameter, namely 
the "ignition" temperature Alternatively we could consider the 

initial laser power density as given and determine the ignition 
temperature required to satisfy the boundary conditions. 

Results for typical conditions represented by a, and m as used by 
Keefer (Ref. 15) are shown in Fig. 6. The initial temperature, T^, was 
{isfiuim'il to ho ')00 K, aiul Lho ignition tempornture. was taken as 9000 K 
representing the onset of ionization for hydrogen at p ■ 1 atm (see Fig. 

5). The temperature distribution through the laser absorption wave is in 
excellent qualitative agreement with results obtained by Keefer (Ref. 15) 
for the same choice of parameters. Of course, since our results represent a 
strictly one-dimensional situation (i.e., radial conduction and finite laser 
beam dimensions are ignored), exact quantitative agreement with Keefer's 
results cannot be expected. 

Another point learned from this one-dimensional exercise is that, with- 
out the explicit inclusion of the radiation loss term in the differential 
equation, the differential equation does not have a solution which is able 
to satisfy the downstream boundary conditions. Both the temperature and its 
derivative will become infinite as x approaches infinity, i.e., the solution 
diverges. Even though the full conservation equations used in the GIM code 
contain the radial conduction term, and therefore provide an energy loss 
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1 atm 
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Fig. 6 - Temperature Distribution Through Absorption Wave (Hydrogen at p - 1 atm) 


loi’m, lliiH iiKiy not bo oL HulLiclunt magtuLudii Lo pruvonL tho toiupurakutro 
Croin reaching unreasonable values. As can be seen ironi the example 
presented by Keefer, the power loss due to plasma radiation far outweighed 
that due to thermal conduction in the radial direction. 

The purpose of investigating this simplified analytical model, although 
it duplicates to some degree work previously accomplished, was to provide 
insights which should be of help to obtain a solution from the GIM code. 

This will be discussed in the following section. 

4.3 TWO-DIMENSIONAL FLUID MECHANICS EFFECTS 

The successful use of detailed numerical methods such as those used in 
the GIM code requires careful planning with regard to grid selection. 

Initial values and boundary conditions, time step selection and possibly 
scaling considerations to accelerate convergence. 

The analytical results shown in the previous section (as well as 
Keefer's results) show that we have to deal with very steep temperature 
gradients in the laser absorption wave. This immediately suggests a rather 
tight spatial grid, but how tight? Selecting too coarse a grid can generate 
oscillations in the solution; but too fine a grid will be costly in terms of 
computation time and storage. 

Using a time-dependent method for the computation of subsonic flow 
fields, particular care must be exercised in order to ensure that boundary 
conditions remain well defined during the course of the calculation. The 
laser heated thruster is modeled as a straight circular duct with vanishing 
velocity at the wall, specified uniform velocity at the Inflow boundary and 
an outflow velocity distribution to be determined by the analysis. In order 
to keep the problem defined, the laser absorption wave must be positioned 
far enough downstream from the subsonic inflow boundary such that the latter 
is not subjected to changes caused by the upstream influence exerted by the 
former. The analytic solution discussed previously shows that, for the 


27 


sample conditions chosen, the laser absorption wave roust be positioned in 
the duct no closer than approximately two cm to the inflow boundary (see 
Fig. 6), assuming that the analytical solution is quantitatively correct 
with respect to gradients expected in the axial direction. 

Initial conditions, and average values for thermodynamic and transport 
properties must be chosen carefully to be representative of the enormous 
temperature range that has to be covered (with constant values) and, simul- 
taneously, they must be consistent in representing the parameters used in 
the calculations with which we want to compare the GIM code results, viz. 



a = 5,0 cm"^ 

2.75 • 10^ cm sec/g 


(33) 

(34) 


A consistent set of variables can be selected as follows; According to 

Keefer's results we expect a maximum temperature of around 20 • 10 K. 

From Fig. 4 it is found that an integrated average value for the specific 

heat covering this temperature range is c” « 110 J/gK. The given ratio of 
— — • _ 2 ^ 

Cp/k then yields k = 4.0 *10 J/cra sec K, a reasonable value as can be 
seen from Fig. 3. While the analytic solutions are determined on the basis 
of global parameters such as pu c^/k, the GIM code requires specification 
of these variables separately. In fact, the GIM code uses y and R to com- 
pute the specific heat. Selecting an integrated average of R » 9.4 J/gK 
from Fig. 5, we can compute the ratio of specific heats to be V* 1.094. The 
equation of state and Eq. (34) then yield 3.593 * 10 g/cm and u » 
50.602 cm/sec, assuming p * 1 atm and T ■ 300 K as initial conditions. 

The selection of the spatial grid for this problem (21 stations in the 
axial direction, 26 points in the radial direction) represents a compromise 
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between anticipated requirenmnta concerning reaolution and the limitations 
of! the MSFG Univac 1100 system with respect to computational speed and data 
storage. 

Pure Clow calculations without laser energy addition were performed 
first. Results from these can be checked against exact (i.e., analytical) 
solutions for axisymmetric Poiseuille flow. Initial calculations with a 
time increment as dictated by the CFL criterion showed extremely slow 
convergence. Much faster convergence was eventually realized by scaling the 
problem. Scaling as applied here involves multiplication of key variables, 
i.e., the velocity, the viscosity and the thermal conductivity in this case, 
by a scale factor (I • 10 , typically) such that Reynolds number and 
Prandtl number are preserved. This type of scaling does increase the Mach 
number. However, because of the rather small velocity, even the scaled Mach 
number remains much smaller than unity, and therefore the basic subsonic 
flow character of the problem is preserved. 

Typical results are showt^ in Figs. ? and 8. As seen in Fig. 7, the 
initial velocity profile is uniform (over the entire spatial domain) except 
for the velocity at the wall which is set to zero. As the parabolic velo- 
city profile developes in time along the channel, its magnitude in the 
outflow plane is continually adjusted by mass balance considerations. This 
adjustment causes a feed-back into the interior of the flow field and 
thereby enhances convergence of the solution. Figure 7 shows the outflow 
velocity profile to have the proper parabolic distribution with a maximum 
velocity on the centerline at roughly twice the value of the average velo- 
city (us represented by the inflow velocity distribution) as predicted by 
the exact solution. A complete velocity vector plot is presented in Fig. 8 
which very clearly shows the contraction and acceleration of the flow in the 
interior of the channel due to the boundary layer developing along the 
channel wall. 

These results also demonstrate well the viscous nature of the flow. 

Note that the scaling factor applied does not affect the balance between 
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vliicouti and inurtiol forcos because the Reynolds number is preserved. Wuen 
Introducing the Healing lactor into the equations it can be seen that its 
oLlect LH to IncroaHe the viscous dissipation tjriiis relative to tiiu convec’" 
tion and conduction terp*,;* in the energy equation. This gives the scaling 
factor basically the function of a convergence accelerator, particularly so 
since it was shown previously that the dissipation terms are negligible for 
very slow flows, and remain so in spite of the scale factor. Typical calcu- 
lations involving a time step of ^t •* I • 10”^ sec, and a thousand time 
steps to obtain a converged solution require about one hour of Univac IlOO 
computer time. 

All efforts to obtain a numerical solution for the flow with laser 
heating were unsuccessful. In order to gain an understanding why solutions 
could not be obtained, it is useful to take a closer look at the time step 
criteria which must be observed to control the stability of the numerical 
calculations. The first one is the CFL condition, which can be expressed as 


^‘^CFL ^ 


u I + a 


( 35 ) 


where A is the smallest spatial grid spacing used, and a is the speed of 
sound. The second is a thermal stability criterion, which, when specialized 
to flows with constant pressure (which is very nearly true for our case), 
reads 



( 36 ) 
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Specializer. to our cate for constant properties and grid spacings of iix * 
0.2 cm and ilr ■ 0.05 cm, these criteria become 



It is immediately evident that At^ presents the more stringent criterion, 
especially for flows with high temperatures as they are expected in the 
laser heated propulsion problem. Since typical cases, using ^t ■ 1 10 ' 

(sec) require about one hour to converge on the Univac 1100, a temperature 
rise to about 20,000 K will reduce the permissible time step by nearly two 
orders of magnitude. Assuming that the total time required to obtain a 
converged solution remains the same as before, we are faced with computer 
times of the order of a hundred hours per case. This is clearly an impos- 
sible proposition. The only remedy here appears to be a computer such as 

I 

the STAR or the CRAY, both of which can better handle the storage require- 
f ments imposed by this problem, and both of which work at much faster compu- 

^ tational speed. Additional computational time savings can most likely be 

I realized by utilizing a variable time step to be evaluated as a function of 

I 

temperature as given by the time step criteria. This feature is included in 
the GIM code version for the STAR computer. Unfortunately, the latter (with 
its GIM code version tailored to it) was not available for this study. 

I 
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5» CONCLUSIONS 


The numerical analysio of two>^dimen»lonal fluid mechanic! effect! in 
the laaer heated thruster la characterised by low denaityt very high temper** 
atures and extremely large gradients in both variables. Because of subsonic 
flow throughout, the governing equations are of the elliptic type and 
require a solution by a time relaxation method with extremely small step 
sizes in time as dictated by the stability criteria. The first prerequisite 
for a successful solution of this problem therefore is a computer facility 
which can handle the storage requirements and is sufficiently fast so that 
the problem can be solved in reasonable time. The Univac IlOO system does 
not meet these requirements. 

Simplified analytical solutions have shown that there exists a partic- 
ular relationship between the gas flow velocity ar;.? the laser power which 
the flow can absorb and simultaneously satisfy the downstream boundary 
condition of a finite temperature. Wlule this relationship is easily 
obtained as part of the analytical solution of the simplified equations, 
most likely it can only be satisfied in a numerical solution via an 
iterative approach. This implies repetitive calculations and therefore 
amplifies the requirements for a large and fast computer. 

Numerical solutions of subsonic flow problems via time dependent 
methods require particular attention to the specification of boundary 
conditions at the inflow and outflow boundaries. This is a subject of 
current research. It is not known at present to what extent the lack of 
appropriate boundary conditions has contributed to the failure of the 
present effort. 

^ The fact that the GIM code version usable on the Univac IlOO system is 
restricted to constant thermodynamic and transport properties should by 
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itaclf not pri>v<^iii(; a aolution. Clearly though, when the tempereture vnriae 
to the oKtent as it does in the Ueer heated thruster problemj the use o£ 
temperature dependent properties would appear to be more realistic* 

Finally, considering the critical dependence o£ the permissible time 
step on the temperature, any future attempts at solving the flow problems in 
the laser heated thruster should use a variable time step methodology in 
order to minimize computer time (and cost). 
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